Analysis of Self-Organized Criticality in the Olami-Feder-Christensen model and in 

real earthquakes 
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We perform a new analysis on the dissipative Olami-Feder-Christensen model on a small world 
topology considering avalanche size difFerences. We show that when criticality appears the Proba- 
bility Density Functions (PDFs) for the avalanche size difFerences at different times have fat tails 
with a q- Gaussian shape. This behaviour does not depend on the time interval adopted and is found 
also when considering energy differences between real earthquakes. Such a result can be analytically 
understood if the sizes (released energies) of the avalanches (earthquakes) have no correlations. Our 
findings support the hypothesis that a self-organized criticality mechanism with long-range interac- 
tions is at the origin of seismic events and indicate that it is not possible to predict the magnitude 
of the next earthquake knowing those of the previous ones. 
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In the last years there has been an intense debate on 
earthquake predictability [1] and a great effort in study- 
ing earthquake triggering and interaction [2-5]. Along 
these lines the possible application of the Self-Organized 
Criticality (SOC) paradigm [6-14] has been discussed. 
Earthquakes trigger dynamic and static stress changes. 
The first acts at short time and spatial scales, involv- 
ing the brittle upper crust, while the second involves 
relaxation processes in the asthenosphere and acts at 
long time and spatial scales [15-21]. In this letter, by 
means of a new analysis, we show that it is possible 
to reproduce statistical features of earthquakes catalogs 
[22,23] within a SOC scenario taking into account long- 
range interactions. We consider the dissipative Olami- 
Feder-Christensen model [12] on a small world topology 
[24,25] and we show that the Probability Density Func- 
tions (PDFs) for the avalanche size differences at different 
times have fat tails with a q-Gaussian shape [26-29] when 
finite-size scaling is present. This behaviour does not de- 
pend on the time interval adopted and is found also when 
considering energy differences between real earthquakes. 
It is possible to explain this result analytically assuming 
the absence of correlations among the sizes (released en- 
ergies) of the avalanches (earthquakes) . This finding does 
not allow to predict the magnitude of the next earthquake 
knowing those of the previous ones. 

The Olami-Feder-Christensen (OFC) model [12] is one 
of the most interesting models displaying Self-Organized 
Criticality. Despite of its simplicity, it exhibits a rich phe- 
nomenology resembling real seismicity, like the presence 
of aftershocks and foreshocks [14]. In its original version 
the OFC model consists of a two-dimensional square lat- 
tice oi N — sites, each one connected to its 4 nearest 
neighbours and carrying a seismogenic force represented 
by a real variable Fi , which initially takes a random value 



in the interval (0, Fth)- In order to mimic a uniform tec- 
tonic loading all the forces are increased simultaneously 
and uniformly, until one of them reaches the threshold 
value Fth and becomes unstable {Fi > Fth)- The driv- 
ing is then stopped and an "earthquake" (or avalanche) 
starts: 
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where "nn" denotes the set of nearest-neighbor sites of 
i. The number of topplings during an avalanche defines 
its size S, while the dissipation level of the dynamics is 
controlled by the parameter a. The model is conserva- 
tive if a = 0.25, while it is dissipative for a < 0.25. In 
the present letter we consider the dissipative version of 
the OFC model with a = 0.21 [30], on a regular lat- 
tice with L = 64 and open boundary conditions (i.e. we 
impose = on the boundary sites). But, in order 
to improve the model in a more realistic way, we intro- 
duce a small fraction of long-range links in the lattice so 
to obtain a small world topology [25]. Just a few long- 
range edges create short-cuts that connect sites which 
otherwise would be much further apart. This kind of 
structure allows the system to synchronize and to show 
both finite-size scaling and universal exponents [24] . The 
curves obtained for different sizes of the system collapse 
into a single one. Furthermore, a small world topology 
is expected to model more accurately earthquakes spa- 
tial correlations, taking into account long-range as well 
as short-range seismic effects [2-5] . In our version of the 
OFC model the links of the lattice are rewired at random 
with a probability p as in the one-dimensional model of 
Ref. [25]. In [24] it was shown that the transition to ob- 
tain small world features and criticality is observed at 
p = 0.02. 
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FIG. 1. (Color online) (a) Power-law distribution of the 
avalanche size S for the OFC model (a = 0.21) on a small 
world topology (open circles) and on a regular lattice 64x64 
(full circles). Data were shifted for clarity. A fitting curve 
with slope T = 1.8 is also reported as full line, (b) PDF of 
the avalanche size differences (returns) x{t) = S{t+ 1) — S{t) 
for the OFC model on a small world topology (critical state, 
open circles) and on a regular lattice (non critical state, full 
circles). Returns are normalized to the standard deviation 
a. The first curve has been fitted with a q-Gaussian (full 
line) with an exponent q ^ 2.0 ± 0.1. A standard Gaussian 
(dashed line) is also reported for comparison. All the curves 
were normalized so to have unitary area. See text for further 
details. 



In Fig.l (a) we plot the distribution of the avalanche 
size time-series S{t) for the OFC model on a small world 
topology (open circles) and on a regular lattice (full cir- 
cles). In our case the time t is a progressive discrete 
index labelling successive events and is analogous to the 
"natural time" successfully used in [21]. We have consid- 
ered up to 10^ avalanches to have a good statistics. In 
both cases the data follow a power-law decay y S^'^ 
with a slope r = 1.8 ± 0.1 even if criticality is present 
only for the small world topology [24]. In the last years 
SOC models have been intensively studied considering 
time intervals between avalanches in the critical regime 
[19]. Here we follow a different approach which reveals 
interesting information on the eventual criticality of the 



model under examination. Inspired by recent studies on 
turbulence and intermittent data [27-29], we focus our 
attention on the "returns" x{t) = S{t + A) — S{t), i.e. 
on the differences between avalanche sizes calculated at 
time t + A and at time t, A being a discrete time interval. 

The resulting signal is extremely intermittent at crit- 
icality, since successive events can have very different 
sizes. On the other hand, if the system is not in a criti- 
cal state this intermittency character is very reduced. In 
Fig.l (b) we plot as open circles the Probability Density 
Function (PDF) of the returns x{t) (with A = 1) ob- 
tained for the critical OFC model on small world topol- 
ogy. The returns are normalized in order to have zero 
mean and unitary variance. The curves reported have 
also unitary area. A behaviour very different from a 
Gaussian shape (plotted as dashed curve) is observed. 
Data are very peaked with fat tails. On the other hand, 
for the model on regular lattice, even if power laws are 
found, the model is not critical since no finite-size scal- 
ing is observed [24]. In this case no fat tails exist, al- 
though a sensible departure from Gaussian behaviour 
is still present (see full circles). These findings sug- 
gest a new powerful way for characterizing the presence 
of criticality. They are also reinforced by similar re- 
sults on other SOC models not reported here for lack 
of space. Another remarkable feature is that such a be- 
haviour does not depend on the interval A considered 
for the avalanche size difference. Also reshuffling the 
data, i.e. changing in a random way the time order of 
the avalanches, no change in the PDFs is observed. The 
data reported in Fig.l (b) for the critical OFC model on 
a small world can be well fitted by a q-Gaussian curve 
fix) = yl[l-(l-g)xV-B]i/(i-«) typical of Tsallis Statis- 
tics [26] . This function generalizes the standard Gaussian 
curve, depending on the parameters A, B and on the ex- 
ponent q. For q — I the normal distribution is obtained 
again, so q ^ 1 indicates a departure from Gaussian 
statistics. The q-Gaussian curve, reported as full line, 
reproduces very well the model behaviour in the criti- 
cal regime, yielding in our case a value of(7 = 2.0±0.1. 
In order to compare these theoretical results with real 
data sets, we repeated the previous analysis for the world 
wide seismic catalog available on line [22] . We considered 
689000 earthquakes in the period 2001-2006. As a fur- 
ther term of comparison, we selected a more complete 
seismic data set, i.e. the Northern California catalog for 
the period 1966-2006 [23]. The latter is a very extensive 
seismic data set on one of the most active and studied 
faults on the Earth, i.e. the San Andreas Fault. In this 
case the total number of earthquakes is almost 400000. 
As pointed out by several authors [14] the energy, and not 
the magnitude, is the quantity which should be consid- 
ered equivalent to the avalanche size in the OFC model. 
In this paper we consider the quantity S = exp{M), M 
being the magnitude of a real earthquake. This quantity 
is simply related to the energy dissipated in an earth- 
quake, being the latter an increasing exponential function 
of the magnitude. 
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FIG. 2. Power-law distribution of 5 = exp(M), M being the magnitude for the World Catalog (a) and for the Northern 
California catalog (c). The correspondent power-law fits are also reported. We show the corresponding values of the magnitude 
M in the upper part of the figures. PDFs of the energy differences x(t) = S{t -I- 1) — S{t) arc shown in (b) for the World Wide 
Seismic Catalog and in (d) for the Northern California Catalog. In both the figures the data have been fitted with a q-Gaussian 
(full line) with q ~ 1.75 ± 0.15. A standard Gaussian is plotted as dotted line in all the figures for comparison. See text for 
further details. 



In Fig.2 (a) and (c) we plot the PDFs of S for the World 

Catalog and the Northern California catalog. Power-law 
decay with exponents r = 2.7 ± 0.2 and r = 3.2 ± 0.2 
reproduce the PDFs for the two cases respectively. Then 
we consider the PDFs of the corresponding returns x{t) = 
S{t + A)- S{t) (with A = 1) and we plot them in Fig.2 
(b) and (d). Also for real data f is a progressive dis- 
crete index labelling successive events. As for the criti- 
cal OFC model previously discussed, fat tails and non- 
Gaussian probability density functions are observed. In 
both cases the experimental points can be fitted by a q- 
Gaussian curve, obtaining an exponent q = 1.75 ±0.15, a 
value which is compatible, within the errors, to that one 
found for the OFC model. The world catalog presents 
large fluctuations in the tails due to the significant in- 
completeness given from the lack of small magnitudes 
events at the global scale. Also for the real earthquakes 
data, by changing the interval A of the energy returns 
X, or by reshuffling the time-series S{t), no change in the 



PDFs is observed. In general both for the OFC and for 

the real earthquakes catalogs the avalanche sizes (ener- 
gies) S occur with a power-law probability p{S) ~ S~'^. 
In the hypothesis of no correlation between the size of 
two events, the probability of obtaining the difference 
x = S'{t + A) - S{t) (whatever A) is given by 

/■OO POO 

P{x) = K dS dS' (5S")"^ 6{S' -S-x) = 
Jo Jo 

dS [S{S+\x\)]-^ (2) 

where K is a normalization factor. The absolute value 
takes into account the exchange of S' with S. The 
integral is divergent for S = 0, so we consider a small 
positive value e as an inferior limit of integration. Then 
one gets 

P{x) = k'-,^—^2F,{t,2t-1;2t;-^-^) , (3) 
zr — 1 e 
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2^1 being the confluent hypergeometric function. The 
probability density function (3) is plotted for various val- 
ues of r in Fig. 3 (a). All these curves can be very well 
reproduced by means of q-Gaussians, whose values of q 
do not depend on e (since we have verified that the lat- 
ter changes only the normalization factor). The relation 
between q and t is shown in Fig. 3 (b), where the points 
are well fitted by a stretched exponential curve. Notice 
that increasing r, i.e. when the power-law tends to an 
exponential, q tends to 1 as expected. The value we get 
for the avalanche size power-law of the OFC model with 
a small world topology is r — 1.8, which corresponds, 
according to Fig. 3, to a value of q ^ 2.1 in agreement 
with the curve shown in Fig. 1(b) within the errors. A 
similar correspondence can be found for the returns of 
real earthquakes data. In particular for r = 2.7 ± 0.2 
and T = 3.2 ± 0.2, see Fig. 2 (a)-(c), one gets values of 
q compatible, inside the errors, with the value q ^ 1.75 
found in Fig. 2 (b)-(d), see the values inside the box of 
Fig.3(b). 
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FIG. 3. (Color online) (a) We show the curves obtained by 
calculating the integral of Eq.(3) considering different values 
of T. (b) We plot (full circles) the relationship between the 
values of q and r, where q is obtained by fitting the curves 
in (a) with a q-Gaussian. The resulting points can be very 

1 —0.795 

well reproduced by a stretched exponential, y = e , 
drawn as dashed line. Inside the plotted box, one can find 
the values of r and q found for the OFC model and the real 
data sets. 



This result explains the q-Gaussian fat-tails in terms 
of differences between uncorrelated (in time) avalanches 
(earthquakes) sizes. On the other hand, we have checked 
that when avalanches are generated by a deterministic 
chaotic dynamics, one finds a dependence on the interval 
considered for the size returns and the resulting PDFs 
cannot be explained with the function (3). 

In conclusion we have presented a new analysis which 
is able to discriminate in a quantitative way real SOC 
dynamics. The results here presented for the OFC 
model and earthquakes data, on one hand give further 
support to the hypothesis that seismicity can be ex- 
plained within a dissipative self-organized criticality sce- 
nario when long-range interactions are considered. On 
the other hand, although temporal and spatial correla- 
tions among avalanches (earthquakes) do surely exist and 
a certain degree of statistical predictability is likely pos- 
sible, they indicate that it is not possible to predict the 
magnitude of seismic events. 
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